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Abstract 

We give a continuum limit value of the lowest moment of a twist-2 opera- 
tor in pion states from non-perturbative lattice calculations. We find that the 
non-perturbatively obtained renormalization group invariant matrix element 
is (x)rgi = 0.179(11), which corresponds to (x)^(2 GeV) = 0.246(15). In 
obtaining the renormalization group invariant matrix element, we have con- 
trolled important systematic errors that appear in typical lattice simulations, 
such as non-perturbative renormalization, finite size effects and effects of a 
non-vanishing lattice spacing. The crucial limitation of our calculation is the 
use of the quenched approximation. Another question that remains not fully 
clarified is the chiral extrapolation of the numerical data. 
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1 Introduction 



Deep inelastic scattering [1] continues to provide important information on the struc- 
ture of hadrons. Phenomenological fits to experimental data give values for the 
moments of parton distribution functions (PDF), including estimates of their er- 
rors, see e.g. [2-4]. Since such moments can be expressed as expectation values of 
local operators, they are accessible to lattice calculations [5]. A direct comparison 
of these lattice calculations of moments with the results of the phenomenological 
fits will test whether these fits are consistent with direct QCD predictions. If we 
think of, e.g., precise determinations of the strong coupling constant from scaling 
violations in deep inelastic scattering, such non-trivial checks are mandatory. 

Lattice results do not come for free, however: Concepts of non-perturbative 
renormalization have to be developed; in the process of moving from the continuum 
of space time to an euclidean lattice a non-vanishing value of a lattice spacing a 
is introduced, leading to discretization effects; running simulations on a computer 
necessitates the use of an only finite volume; the limited amount of computing 
resources leads to the fact that simulations are performed at rather large values of 
the quark masses that are far from their values assumed in nature. In addition, 
at present we are left with the quenched approximation, neglecting internal quark 
degrees of freedom. Finally, the numerical results are plagued by statistical errors 
that can be substantial for bad choices of operators. 

In the course of our work [6-11] to reach a reliable value for a moment of a 
parton distribution function we eliminated important sources of systematic errors 
besides the quenched approximation (see [12-14] for summaries of these works). The 
transition to full dynamical simulations is undertaken world-wide today [15] and the 
next years will see the exciting results of such calculations. Another open question is 
the chiral extrapolation that is not understood presently (see [16] and refs. therein). 

Let us sketch how we have addressed the systematic errors of the lattice calcu- 
lations in our work: 

• Non-perturbative renormalization 

We adopt in our work [6,7,11] the Schrodinger functional (SF) scheme [17], 
which has been proved to be a very successful and practical way to compute 
scale dependent renormalization constants. By evolving deep into the pertur- 
bative regime, renormalization group invariant (RGI) quantities can be deter- 
mined that allow to relate lattice results into any desired continuum renormal- 
ization scheme [18-20]. In this paper we compute the renormalization factor 
at the matching scale (see sect. 13.21 for details). 

• Discretization effects 

We have controlled effects of a non-vanishing lattice spacing a by performing 
all our lattice calculations with several values of a with two different lattice 
formulations of lattice QCD: one is the standard Wilson fermion, the other 
the non-perturbatively improved Wilson fermion formulation. In this way, all 
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quantities could be extrapolated to their continuum value in a controlled way 
performing a constrained contimium limit [8,11]. In particular in this work 
we have used this strategy to control the continuum limit of the renormalized 
matrix element. 

• Finite Size Effects 

The finite size effects were controlled in two ways. For the computation of 
the evolution of the renormalization constants [7, 9, 11], the finite volume 
Schrodinger functional scheme was used, in which the scale is identified 
with the inverse lattice extent L"^ . In this way, the finiteness of the lattice 
has been utilized to determine the scale dependence. 

For the matrix elements, a careful test of finite size effects has been performed 
[21] with the somewhat surprising result that such effects are rather large for 
the pion matrix clement at a point where the pion mass itself shows no effects. 
Since this was done for pion matrix elements only, it will be very important to 
repeat the analysis for baryon matrix elements, where such effects are expected 
to be even larger [22]. 

• Statistical errors 

By employing generalized boundary conditions in space, the signal to noise 
ratio of operators that need an external momentum could be minimized by 
optimizing parameters that arc introduced by the generalized boundary con- 
ditions [11]. This allowed us to obtain results that would not have been possible 
otherwise. 

• Chiral extrapolation 

The extrapolation of the numerical data from the rather heavy quark masses, 
where the simulations are performed, to their physical values is still under 
debate and not clarified, see [16] for a general discussion, and [23-28] for 
works on the chiral extrapolation of PDF. We will follow the strategy in this 
work that we first extrapolate the values of our non-pcrturbatively improved 
matrix elements to the continuum limit. Comparisons to predictions of chiral 
perturbation theory are then performed directly in the continuum in order to 
avoid lattice spacing effects that would render the interpretation difficult. 

After having understood and overcome the above difficulties, we can now pro- 
vide results that are free of these systematic uncertainties besides the, presently still 
unavoidable, quenched approximation. In this work we summarize our result and, 
most importantly, we provide the missing part, i.e. the renormalized matrix element 
at the hadronic scale /iq — 275 MeV down to which the scale evolution has been 
computed ( [11]), with a full control on the continuum limit and a carcfull analysis of 
the chiral extrapolation. As the most important quantity we consider the renormal- 
ization group invariant (RGI) matrix element of the twist-2 operator between pion 
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states corresponding to the first moment of tlie valence quark distribution, wliicli 
we find 

(x)rgi = 0.179(11). (1) 

This quantity is of central importance since it allows to relate the non-perturbati- 
vely obtained results using a particular lattice renormalization scheme to more con- 
ventional schemes. For example, if we use the MS scheme, we find that at a scale of 
2 GeV the value of the matrix element compares to phenomenological estimates of 
the same quantities extracted from global fits of experimental data [29] as follows *: 

{x)^{fi = 2 GeV) = 0.246(15) 
(xW(/i = 2GeV) = 0.21(2). (2) 

This results demonstrates that the running of moments of parton distribution 
functions in the continuum can be calculated from lattice simulations in an intrinsic 
perturbative scheme like the MS retaining all the non-perturbative information com- 
ing from the RGI matrix element. In addition, the errors coming from the lattice 
simulations are comparable to the experimental errors, which opens the possibility 
to perform direct tests of QCD as the theory of the strong interaction. 

2 Renormalized matrix element 

The moments of parton density distributions are related to expectation values of 
local operators, which are renormalized multiplicatively by applying appropriate 
renormalization factors Z{^) that depend on the energy scale fj, (see, e.g. ref. [1]). 
This leads to consider renormalized matrix elements 0'^™(/i) using some renormal- 
ization scheme, denoted by ren. For lattice calculations, we are aiming at here, 
a very useful scheme is the Schrodinger function (SF) renormalization scheme [17] 
since it applies in small volumes. 

If the energy scale fj, of Z{fi) is chosen large enough, it is to be expected, and 
indeed it can be checked explicitly, that the scale evolution is very well described 
by perturbation theory, giving rise to the following definition of a renormalization 
group invariant (RGI) matrix element: 

Orgi = 0'^i^^) ■ f'ifif^)) (3) 

with g{^) the running coupling computed in the same SF scheme and 

*Some details on the extraction of the phenomenological number can be found in sect. IS.ll of 
this paper. 
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where l3{g) and '-f{g) are the (3 and anomalous-dimension functions computed to a 
given order of perturbation theory in a specified scheme, i.e. here the SF scheme. 
Once we know the value of Orgi evaluated non-perturbatively, the running matrix 
element in a preferred scheme can be computed, for example in the MS scheme: 

o^{^^) = o^oi/f^igU^^)) (5) 

with now, of course, the (3 and 7 functions computed in the MS scheme. Thus, 
although the SF is an unphysical finite volume scheme (but therefore most suited 
for lattice simulations) it can be related to more conventional continuum renormal- 
ization schemes by providing the renormalization group invariant quantities. 

A non-perturbatively obtained value of the renormalization group invariant ma- 
trix element is hence of central importance. Its calculation has to be performed in 
several steps. The reason is that we have to cover a broad range of energy scales 
- from the deep perturbative to the non-perturbative region. Using the scale de- 
pendent renormalization factor we write the renormalized matrix element 
of eq. as 

where \h) is the hadron state we are interested in. So far, all our discussions have 
been in the continuum. However, if we think of the lattice regularization and even- 
tual numerical simulations to obtain non-perturbative results, it would be convenient 
to compute the renormalized matrix element at only one (i.e. small hadronic) scale 
/iQ. We therefore rewrite the r.h.s. of eq. (jH)) as 

{h\0\h) {h\0\h) Z^^i^^o) 

where we introduce the step scaling function a{fi/ fio,g{fi)), which describes the 
evolution of the renormalization factor from a scale fiQ to a scale /i. The advantage 
of concentrating on the step scaling function instead of the renormalization factor 
itself is that the step scaling function is well defined in the continuum and hence 
suitable for eventual continuum extrapolations of lattice results. We finally write 
the r.h.s. of eq. as 

Orgi = OS^(/io) a(/^/^o, gif^)) ■ f'ififi)) (8) 

v ' 

= ®INV (mo) 

with 0^^(/io) the renormalized matrix element, which is to be computed only once at 
a scale /xq and the ultraviolet (UV) invariant step scaling function (3jJ^^^(/io), which 
still depends on the infrared scale fio, and on the renormalization scheme adopted. 
In ref. [11] we have given a value for the UV invariant step scaling function and 
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checked the independence on the ultraviolet scale In this work we will provide 
the missing part, i.e. the rcnormalized matrix element at the scale /xq- The final 
result for the renormalized matrix element of the pion that we are going to provide 
here, was made possible by a number of theoretical and conceptual developments 
that we could achieve over the last years [6-9]. The methods developed there can 
immediately be taken over to other matrix elements than the lowest twist case of 
the pion considered here and to the unquenched situation. Nevertheless, with this 
paper we want to finish the analysis of the pion matrix element with the aim to have 
eliminated important systematic uncertainties besides the quenched approximation. 



2.1 Transfer matrix decomposition 

The moments of parton distribution functions (PDF) are related to matrix elements 
of leading twist r (r = dim - spin) operators of given spin, between hadron states 
h{p) 

{h{p)\0,,...,Jh{p)) = M(^-i)(^)p,,---p^„ 

+terms 5/,,^., (9) 



(pi) = m(^-^) (a* = g) . (10) 

We concentrate in this work on the twist-2 operator corresponding to the second 
moment of the parton distribution functions (PDF) between charged pion states. 
In the following we consider the fermionic fields ip{x) and ijj{x) as doublets in the 
flavour space. In particular we will concentrate on the valence u or d distribution 
as explained below. This amounts to consider operators of the form 

1 - ^ 

Ofj,u{x) = -^{x)-f{^ Du} i^ix) - 5^^ ■ trace terms , (11) 
where {• • • } means symmetrization on the Lorentz indices, 0^=0^ — and 

Bm=^(V^ + V;), 5,= ^(V,+ V*). (12) 

The definitions of the lattice derivatives and conventions are given in the appendix. 

There are two representations of such an operator on the lattice [30]. The first 
representation takes /x 7^ i^, whereas the second uses n = u. The precise definitions 
of the operators used here are 

Ci2(x) = -Mx)^{i B2} il^{x) (13) 



and 



1 - r - 1 ^ - 1 

044(2:) = -^{x) 74 Da -o 5^7fe Dk i^{x) . (14) 
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In computing the matrix elements of these operators, a non-zero momentum in two 
different spatial directions has to be supplied for Oi2{x) in eq. (fT^ . whereas for the 
operator 044(0;) in eq. (fT^ no momentum is needed. It is to be expected, and indeed 
verified in numerical simulations, that the signal of the matrix element of 044 (x) is 
better than the one of Ouix). Thus in the following investigation we consider only 

Our setup of lattice QCD is on a hyper-cubic euclidean lattice with spacing a and 
size X T. We impose periodic boundary conditions in the spatial directions and 
Dirichlet boundary conditions in time, as they are used to formulate the Schrodinger 
functional (SF) [17] (we refer to these references for unexplained notations). Using 
homogeneous boundary conditions, where the spatial components of the gauge po- 
tentials at the boundaries and also the fermion boundary fields are set to zero, the 
Schrodinger functional partition function can be written as [31] 

Z = {^o\e-™F\^o) , (15) 

where the initial and final states \io) carry the quantum numbers of the vacuum 
and P denotes a projector on gauge invariant states. The states with charged pion 
quantum numbers in the Schrodinger functional are, indicating with ( and ( (and 
the corresponding (' and (') a flavour doublet, the dimensionless fields 

^ = |l E C(y)75r+C(z) and S' = ^ C'(u)75r-C'(v) , (16) 

y,z u,v 

where = -^(r^ ± ir^) and r'^ with k = 1, 2, 3 are the usual Pauli matrices. The 
pion interpolating fields § and §' are respectively localized at a;o = and xq = T. 
The desired matrix element is obtained from the correlation function 

/44(Xo) = -i(§'044(x)§), (17) 

where we have used the independence on the spatial components x of x = (a;o,x). 
The Wick contractions of this correlation function contain also a disconnected piece 
that we neglect consistently with the fact that we are interested on valence quark 
distribution. For normalization purposes it is important to define the boundary to 
boundary correlation function 

/i = -^(S'S)- (18) 

The basic fermionic Wick contractions for these two correlation functions are de- 
picted in fig. m 

In the following we will denote with \io) and the states carrying respectively 
the quantum mumbers of the vacuum and of the charged pion at zero momentum. 
The correlation functions /44 and /i have the following quantum mechanical repre- 
sentations 

/44(xo) = Z-il(v|e-(^-^°)=P 044 e-^°=P|^.) , (19) 
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Figure 1: Fermionic Wick contractions for the /44 and /i correlators. C and C denote the 
boundary gauge field respectively at a;o = and a;o = T; x denotes the insertion point of the local 
operator. 



h = Z-']^{i^\e-™¥\i^) . (20) 

In order to extract the pion mass we have analyzed also the improved axial correla- 
tion function 

/iK) = -yK(^)S), (21) 

where 

Al,{x) = Ao{x) + acA\{dl + d^)P{x) ; (22) 

ca is the improvement coefficient and the axial and pseudoscalar local operators take 
the form 

Aq{x) = '0(a;)7o75T"'0(a;) (23) 

P{x)^i^{x)^^T-^l^{x). (24) 

The non-perturbative value of ca was taken from ref. [32]. 

Following ref. [31] we insert a complete set of eigenstatcs of the hamiltonian 
and, retaining only the first non-leading corrections, we have for the improved axial 
current correlation function /i(a;o) and /i, 

r3 

/j^(xo) - — p(0, 0|A^|0, 7r)e-"*-^o{l + Ty^e"^"^ + Ty^e-^^-^o)'"^} (25) 



2 



h ^p^e-''^ . (26) 
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where the coefficient appearing in this expansions are 



P 



(0,7r|v) 

(0,0|2o)' 


(27) 


(0,0|A^|l,7r)(l,7r|v) 


(28) 


(0,0|A^|0,7r)(0,7r|z^) ' 


(zo|l,0)(l,0|Ao|0,7r) 


(29) 


(2o|0,0)(0,0|Ao|0,7r) ■ 



The energy gap in the pion channel between the fundamental and the ffist excited 
state is denoted by A and is estimated to be Aro ~ 3.2, while is the mass of the 
0"*""'" glueball, mcro ~ 4.3 [31]. For the matrix element we find 

/44(xo) ^ ^p'(0,vr|O44|0,vr)e-'"-^{l + r^5^^e-^'°^ + r^S44e"^^"^°^^} , (30) 
where we define the ratio 



(0,7r| 


1^44 1 


|l,7r)(l,7r| 


\Q 


(0,vr| 


1^44 1 


|0,7r)(0,7r| 


\in) 



^5.4= r:;.^,;,;;; . (31) 



A corresponding transfer matrix decomposition is obtained for the correlation 
function /12. From these expressions it becomes clear that the matrix element we 
are interested in, neglecting contributions from excited states, can then be extracted 
from the plateau value of the following ratio: 

^^ = {0,n\O,,\0,rf). (32) 
Ji 

Finally, in order to relate this numerically computed ratio with the corresponding 
continuum operators in Minkowski space, we need a suitable normalization factor 

M = — (0,7r|C»44|0,7r) (33) 

with K the standard hopping parameter of the Wilson-Dirac-operator. We remark 
here that {x) corresponds to the valence distribution of a single quark {u or d for 
example). We followed ref. [31] in order to extract the plateau values for the effective 
pion mass and the matrix elements. Using the transfer matrix decomposition in 
eq. ()25|) . the effects of higher excited states for the effective mass are given by 

me//(xo) ~ m^+A t]^ e"^ ^"-mc Va e"^" ^^""""^ • (34) 

In fig. 121 we show the effective pion mass as function of the anticipated excited state 
contaminations given in eq. ()34|) . From the linear behaviour of the effective mass as 
a function of ^° and e"™''^ (t-xq) conclude that these are indeed the leading 
corrections. 
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Figure 2: Effective pion mass as a function of the theoretically expected excited state contribu- 
tions obtained from a transfer matrix decomposition. Data corresponds to the simulation point at 
l3 — 6.45 for three k values with the non-perturbatively improved action (see table In the left 
panel we plot the effect of an excited pion state, whereas in the right panel we show the effects of 
the glueball contribution. 

The contribution from excited states on the plateau value of the matrix element 
following eqs. (|^ and (plj) is given by 

(x)(xo) ^ (x) + e"^ (^--))} . (35) 

We show in fig. El the matrix element as a function of the expected excitation in 
eq. ()35p . Again we observe a linear behaviour of the matrix element indicating 
that the corrections come mainly from the first excited pion state. It is important 
to remark that itlg and A have been computed using different boundary conditions 
(see ref. [31]), where in principle excited states corrections have different amplitudes. 
The agreement between the data and the expected form is then reassuring that the 
excited states contamination is well controlled. 

From fig. El and fig. 01 we can read off, what is the systematic error on the pion 
mass or the matrix element value if a particular value of the time separation xq is 
taken to extract their values. In our analysis, we demanded a systematic relative 
error, for all the (3 and k values, coming from the excited states of 0.1% for the pion 
mass and 0.4% for the matrix element which is well below the statistical accuracy 
of our computations. Relating the value of Xq to physical units, the above choice for 
our desired accuracy leads to a window for the extraction of the plateau that for all 
the f3 and k values are around 

- 1.2 fm < a;o ^ T-1.1 fm for rUeg, 
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Figure 3: Matrix element as a function of the excited state contribution as obtained from a 
transfer matrix decomposition. Data corresponds to the simulation point at /3 = 6.45 for k = 0.1350 
with the non-perturbative improved action (see table EJ- In the plot we superimpose the results 
for XQ > T/2 and xq < T/2. 

- 1.3 fm < Xq ^ T-1.3 fm for the matrix element. 

We give in fig. 0]an example for the plateau behaviour of the effective mass and in 
fig. Elan example for the plateau behaviour of the matrix element. In tabled we 
summarize the time intervals chosen for m,r and (x) that fulfill the aforementioned 
conditions. 

3 Numerical details and results 

In this section we will give numerical details and results about the computation of 
the bare matrix element and of the renormalization constant at the low energy scale 

fiQ. 

3.1 Bare matrix element 

We have performed a set of quenched simulations for five j3 values, varying the lattice 
spacing between a = 0.093 fm and a = 0.048 fm, see table [H 

In order to have a better control over the continuum extrapolation of our lat- 
tice results we performed two independent sets of simulation at these /? values, 
one employing standard Wilson fermions and the other using non-perturbatively 
0(a)-improved Wilson fermions [32,33]. As mentioned already above, we used the 
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Figure 4: Plateatix for the effective pion mass. The fit region to extract the mass is indicated 
as a solid line. The simulation was done at /3 = 6.45. 
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Figure 5: Plateau for the matrix element. The fit region to extract the mass is indicated as a 
solid line. The simulation was done at /3 = 6.45. 
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Lja 


T/a 


N 


6.0 


5.37(2) 


16 


32 


600 


6.1 


6.32(3) 


24 


42 


600 (Wilson) 


6.1 


6.32(3) 


32 


56 


185 (Clover) 


6.2 


7.36(3) 


24 


48 


600 


6.3 


8.49(4) 


24 


64 


400 


6.45 


10.46(5) 


32 


72 


400 



Table 1: Parameters of our simulation points. N denotes the number of measurements taken 
into account. For all /3-values except for /3 = 6.1 the number of measurements refers to Wilson 
and Clover simulations individually. 

quenched approximation throughout this work. In order to achieve an extrapolation 
to the chiral limit, we employed three values of the quark mass corresponding to 
a value of that lies in the range of 550 MeV - 1 GeV for all the /3 values. For 
the simulations of the lightest quark mass at each /5 value (the corresponding pion 
masses range in 4.5 < vtlt^L < 5.2), we have corrected for the finite size effects, 
following ref. [21]. For all the /3 values the finite size corrections are in the region 
0.5% — 1.3%. The summary of our bare results in tableElare thus free from finite size 
effects. The systematics coming from the higher excited states have been controlled 
as explained in the previous section. 

3.2 Renormalization constant 

In order to renormalize the bare matrix element at the scale /iq, where we can 
make contact to the running described by the non-perturbatively computed UV 
invariant step scaling function ©|^^^(/io), we have to compute the renormalization 
constants Z^^(/io). The continuum limit of the renormalized matrix element requires 
to compute Z^^(yUo) at exactly the lattice spacing, where the matrix element has 
been calculated, while keeping the scale fixed. Decreasing the lattice spacing a, we 
would hence have to increase the lattice volume in order to stay at a fixed value of 
IIq = (1.436ro)~^. Since we can not vary the lattice size continuously, we performed 
instead simulations on a sequence of lattice sizes and adjusted the values of /? to 
realize the correct value of no- The values of (3 are slightly different from the ones 
used in ref. [34], and were obtained [35] using a new determination of r^/a [36]. 
We recall that for our determination of the renormalization constant at the scale fiQ 
using the finite volume SF scheme we have performed simulations directly at k = Kc, 
determined through the PCAC Ward identity. Details about the computation of the 
renormalization factor, such as renormalization condition and external parameters 
typical of the SF scheme, can be found in ref. [11]. We give in tableElthe parameters 
of our runs and the values of the Z factors at fiQ for both Wilson and 0(a)-improved 
Wilson fermions. 

In order to match the (3 values where the matrix elements themselves have been 
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Figure 6: Numerical results for Z^^^{gQ , hq) together with their interpolating polynomials, 
computed, we describe the P dependence by an effective interpolation formula, 

2 

Z'^if,o)=J2zf{(3~6.0y (36) 

1=0 

with the coefficients for Wilson and 0(a)-improved Wilson fermions listed in table 
131 The statistical uncertainty to be taken into account when using this formula 
is about 1.1% for the non-perturbatively improved clover action and 1.4% for the 
Wilson action ^ . In fig. IHlwe show the P dependence of Z^^{^q) for the operator 
and for the two actions used, together with the plot of the interpolating formula 

Using this interpolation formula, we can match the values of a corresponding 
to the (3 values used for the computation of the bare matrix element, which allows 
finally to obtain the renormalized matrix element at the scale /Xq in the continuum 
limit. 

We have performed a continuum extrapolation at fixed values of {rQiTiT^)'^. Fixing 
(romTr)^ is achieved by interpolating our bare matrix elements linearly in the quark 
mass. The physical values of the pion mass range 550 MeV < < 1 GeV. In 
fig. m we show an example of the continuum limit at our next to lowest pion mass 

■fFor (3 « 6.5 the error to be associated grows to 1.4% for the non-perturbatively improved 
clover action and to 2.0% for the Wilson action. 
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Figure 7: Combined continuum extrapolation for the pion matrix element computed with Wilson 
and 0(a)-improved Wilson fermions, at (rom^)^ = 2.57, which corresponds to a pion mass of 
= 632 MeV. 



((rom^)^ = 2.57). The continuum extrapolation of the renormahzed matrix element 



lim 

a->0 



{n\OM\n){a,rom^) 



(37) 



has been done via a constrained linear fit in the lattice spacing a using simulation 
results obtained with the Wilson and clover action. For all the values of the pion 
mass we have performed this constrained extrapolation excluding the coarsest lattice 
of our data set. In sectional we will discuss the uncertainties connected with the 
continuum extrapolation, and the one related to the chiral extrapolation. If we 
perform the chiral extrapolation linearly in the pion mass squared directly in the 
continuum we obtain a value of 



(x)^^(/io) = 0.810(33) 



(3^ 



for the continuum renormalized pion matrix element in the SF scheme. 

We can anticipate here for clarity that, in this particular case, to invert the order 
of the chiral and continuum extrapolation gives a fully consistent result, indepen- 
dently on how the continuum extrapolation is performed. 
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4 RGI matrix element 



The phenomenological analysis of the experimental data is usually done in the MS 
scheme. In order to translate our fully non-perturbatively computed matrix element 
(x)'^^(/io) renormalized in the SF scheme to the MS scheme, we need first to calculate 
the universal renormalization group invariant matrix element. This is done following 
the complete non-perturbative evolution [11] of the step scaling function in the 
SF scheme as outlined above. Using the UV invariant step scaling function as 
determined in our earlier work [11], it is now possible to eliminate any reference 
to the scale /iq and to obtain finally the RGI matrix element. In ref. [11] we have 
already computed 

©7w^(/io) = 0.221(9) . (39) 
Thus using formula (jH)) we can obtain the main result of this paper 

{^)rgi = i^f^if^o) Q^Nfil^o) = 0.179(11) . (40) 

The RGI matrix element allows a simple conversion to any desired scheme 
(e.g. MS at 2 GeV) requiring only the knowledge of the /3- and 7-function up to 
a certain order in perturbation theory in this scheme [37,38]. The total renormal- 
ization factor to directly translate any bare matrix element of into the RGI 
matrix element, can be written as 

^RGi(^7o) = ^uv,SFt . ■ (41) 

^INV l/^oj 

Combining eq. ()39p and the interpolating formula (jHUj) we obtain a further interpo- 
lating polynomial 

2 

^RGi = 5^^r^(/5-6.0)^ (42) 

whose coefficients are listed in table 01 These parametrizations of ^rgi are to be 
used with the same uncertainty of Z^^ and an additional error of 4.0%, coming from 
the uncertainty of &^^y^ (hq), that has to be added quadratically after performing a 
continuum extrapolation. Our resulting number (x)'^^(/i=2GeV)=0. 246(15) is not 
fully compatible with previous lattice computation [40], where it was found that 
(x)^^(/i=2.4GeV)=0. 273(12), and in better agreement with the result coming from 
the global fits [29] of experimental data {x)phen{f^ = 2GeV) = 0.21(2). The disagree- 
ment of these numbers could be explained by the fact that in this paper we analyze 
and correct for several sources of systematic errors. We apply a non-perturbative 
renormalization and perform the continuum limit, while in [40] a perturbative renor- 
malization was adopted without performing the continuum limit (only one lattice 
spacing). Moreover we correct, where present, for finite size effects. 
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5 Chiral and continuum extrapolation 



Since the form of the correct chiral extrapolation of the lattice data obtained at 
rather large quark masses is still under debate (see [16] for a recent review and refs. 
therein, while for the problem considered here see [23-28] ), we aimed at avoiding 
a possible source of systematic error and first performed a continuum extrapolation 
and then tried to extrapolate the data to the chiral limit. In this way, continuum 
chiral perturbation theory is applicable and we do not have to worry that chiral 
symmetry breaking effects may spoil the chiral extrapolation at non- vanishing values 
of the lattice spacing. We fixed the physical values of the pion mass in units of tq [39] 
and performed a continuum extrapolation keeping rom.^ fixed. To this end we had 
to slightly interpolate the values of the matrix elements, since in our simulations 
the pion mass at different values of f3 were not obtained at the same value of the 
physical pion mass. The corresponding error, taking the correlations into account, 
of this slight interpolation is, however, well below the statistical accuracy of our 
numerical data. We have performed the continuum limit of the matrix element from 
our simulation points using three, four and five values of /?, employing Wilson and 
0(a)-improved Wilson fermions, in order to study the effect of including coarser 
lattices for the continuum extrapolation. This detailed analysis is summarized in 
fig. IHl In all the following discussion, when we talk about chiral extrapolation, we 
use a linear extrapolation in the squared pion mass. 

In this figure we plot the value of {x)^^{fi = 2 GeV) as a function of the number 
of points used to perform the continuum extrapolation. Moreover, for each choice 
of the number of points used, we compare several ways to perform the continuum 
and chiral limit (the symbols refer to the symbols in fig. |S)): 

1) first the constrained continuum limit and then the chiral extrapolation (V); 

2) first the chiral extrapolation and then the constrained continuum limit (A); 

3) first the chiral extrapolation, then the cont. limit using only the Wilson data (□); 

4) first the chiral extrapolation, then the cont. limit using only the clover data (O). 
All these lattice results are compared with the same quantity obtained through 

global fits of the experimental data [29]. Our conclusions, for this quantity and 
for the range of masses and lattice spacings simulated are twofold. First of all, 
the order of how to perform the chiral and continuum extrapolation is irrelevant, 
independently from which action is used to do the continuum limit. Then the 
continuum limit is consistent if performed using three, four or five points, also in 
this case, independently of the lattice action used. 

We decided to perform first the continuum limit at fixed pion mass using the 
four points corresponding to the four smallest lattice spacings. Only as a second 
step we then perform the chiral extrapolation linearly in the squared pion mass. 
The continuum extrapolated data at fixed are collected in table El 

A remark is in order here. It has been shown in ref. [25] that in quenched chiral 
perturbation theory (Q^PT) at the leading order, this particular matrix element 
is free from chiral logarithms. In ref. [27] it has been shown that in the full xPT 
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Figure 8: Value of {x)^^^{iJ. = 2 GeV), obtained from our lattice simulations, as a function of 
the number of points used to perform the continuum extrapolation. For each choice of the number 
of points used we compare several ways to perform the continuum and chiral limit: (V) first the 
constrained continuum limit and then the chiral extrapolation; (A) first the chiral extrapolation 
and then the constrained continuum limit: (□) first the chiral extrapolation and then the continuum 
limit using only the Wilson data; (O) first the chiral extrapolation and then the continuum limit 
using only the clover data. The band represents the same quantity, with its error, obtained through 
global fits of the experimental data [29] . 
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Figure 9: Chiral extrapolation in the continuum. The dashed curves refer to the formula 143|l 
with A = {0.4,0.7,1.0} GeV and = 93 MeV. The chiral extrapolated value (V) is compared 
with the phenomenological estimate (^). 



indeed there are chiral logarithms. In ref. [24] an effective phenomenogical ansatz is 
given by 



1 2^ ^ 



In 



Bnmi , (43) 



motivated by the chiral expansion in the full theory. A similar formula was suggested 
to explain the discrepancy between the lattice results and the global fits of the 
experimental data of the proton (x). 

In fig. El we show the continuum values of the renormalized matrix elements as a 
function of the pion mass (see table El). In order to illustrate possible unquenching 
effects we also perform a chiral extrapolation using eq. ()43|) . For the different dashed 
curves in fig. Elwe used different values of the parameter A = {0.4, 0.7, 1.0} GeV that 
is supposed to provide an estimate of the size of the pion cloud. However, the data, 
even in the continuum limit, clearly prefer a straight line extrapolation consistent 
with the results of QxPT. Using a linear chiral extrapolation we obtain a final value 
for our matrix element consistent with the experimental results, up to quenching 
effects. In order to disentangle the quenched and chiral uncertainties, lattice data 
have to be provided very close to the physical point, where the pion assumes a mass 
as measured in experiment; presumably chiral invariant formulations of QCD are 
necessary to clarify the question of the chiral extrapolation. For the time being and 
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for this work we take the straight hne behaviour as suggested by the numerically 
obtained data for the chiral extrapolation. 

5.1 Discussion and concluding remarks 

The main results of our paper, i.e. the renormalization group invariant matrix ele- 
ment and the renormalized matrix element at a scale of 2 GeV in the MS scheme, 
are summarized in eqs. ()1I2|) . Comparing our results to phenomenological estimates, 
it comes out to be slightly higher, but already more precise than the latest NLO 
analyses of Drell-Yan and prompt photon ttN data [29] , which yield a combined ex- 
perimental value of (x)'^^(/i=2GeV) = 0.21(2). These results have been computed 
for the valence quark distribution ^. We conclude that the non-perturbatively ob- 
tained number for the twist-2 matrix element of a pion is not completely consistent 
with phenomenological estimates. Although this should not be too worrisome, since 
we are still left with the quenched approximation, it is worthwhile to discuss the 
lattice numbers further. 

One cause for the deviation could lie in the chiral extrapolation of this matrix 
element that was done linearly in the pion mass. If we would use the phenomeno- 
logical fit ansatz of ref. [24], we find (x)^(2GeV) = 0.21(1)(2), where the first 
error corresponds to varying the data in their errorbars and the second error is the 
variation for A in the range 0.4 to 1.0 GeV. Clearly this chiral extrapolation would 
reconcile the lattice results with experiments. However, as can be seen in fig. El at 
the moment it is not possible to test such an ansatz unambiguously, since the values 
of quark masses used in the simulations are too large. Only with future simula- 
tions, presumably employing e.g. chiral invariant formulations of lattice QCD, this 
question could be clarified. 

Another cause of the deviation could certainly be the quenched approximation, 
which is a severe limitation and effects of dynamical quarks have to be explored 
in the future. In any case, we conclude that the Schrodinger functional method 
has been proven to provide an excellent tool to compute the non-perturbative scale 
evolution of multiplicatively renormalized quantities such as the quark mass [34] or 
the matrix elements [11] considered in this work. 

In order to further utilize our results, we attempted to use the non-perturbatively 
obtained renormalization constants to see the effect of a non-perturbative renormal- 
ization on the results of other groups on the quenched bare matrix element of the 
O44 operator between proton states (prehminary results can be found in [41]). 

The LHPC coll. [43] has computed the bare matrix element of the operator 
between proton states at /5 = 6.0 using the Wilson action. The QCDSF coll. [42] 
has computed the bare matrix element of the 044 operator between proton states 
at /9 = 6.0, 6.2, 6.4 using the non-perturbatively improved clover action. Both the 
collaborations perform a chiral extrapolation linear in the quark mass at fixed lattice 

■I- We remark here that, as stated in ref. [29], the extraction of the sea distribution in a pion is 
currently impossible due to the lack of the necessary data at low Bjorken x. 
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Figure 10: Continuum extrapolation of the non-perturbatively renormalized matrix element of 
O44 between proton states based on the quenched bare data, both already chirally extrapolated 
linearly in the quark mass, from refs. [42] and [43], using the non-perturbatively computed Z factor 
of the present work (ZeRo coll.) 
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spacing. If we take now the non-perturbative renormalization factor computed in 
this paper (ZeRo colL), we obtain the resuh summarized in fig. where we show 
the continuum extrapolation of the non-perturbatively renormahzed matrix element 
of C44 between proton states based on the bare data, already chirally extrapolated, 
from refs. [42], and as a comparison the same non-perturbatively renormahzed ma- 
trix element using the bare data of ref . [43] . In fig. we also indicate the result 
obtained from global fits of the experimental data as listed in ref. [43] (see refs. 



inside for details). The final numerical results are 

{x)^{l2 = 2GeV) = 0.227(14) ZeRo + QCDSF (44) 

(x)^(/i = 2GeV) = 0.243(20) ZeRo + LHPC /3 = 6.0 (45) 

(x)^(/i = 2GeV) = 0.154(3) global fits (46) 



For the first time a continuum extrapolation in the case of the proton matrix element 
could be done, thanks to the non-perturbative renormalization factor computed in 
this paper. 
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Appendix 

For the computation of the matrix elements the quark fields are chosen to be periodic 
in the three space directions, 

ipix + Lk) = ip{x), %l){x + Lk) =%l){x) . (47) 

The lattice derivatives in the forward direction are given by 

V^^/'(x) = '^[U{x, /i)^(a; + afi) - ^{x)] (48) 

V*?/'(x) = -[i^i^x) - U{x - afi, fi)~'^tlj{x - afi)] (49) 
and the backward derivatives are defined by 

^(x) = -[^(x + afi)U{x, /i)-^ - i^{x)] (50) 
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a 



ip{x — ajj,)U{x — aji, fj,)] . 



(51) 
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Tabic 2: Rcs^ilts for the pscudoscalar mass and the bare matrix element at aU our simulation 
points. We specify also the fit interval for the effective mass and the bare matrix element obtained 
from the request of having the systematic uncertainty coming from the excited states below 0.1% 
and 0.4%, respectively. See text for further details. 
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Table 3: Results for Zq^^ and Z^^^ at the scale /^o = (l-436ro) ^ with Wilson and non- 
perturbatively improved clover actions. 
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Table 4: Coefficients of the interpolating polynomials, eq. (|36|l and (|42(l . Uncertainties are 
discussed in the text. 
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Table 5: Results for the matrix element in the continuum limit (quenched approximation), con- 
verted to the MS scheme at 2 GeV. The continuum extrapolation was performed with a constrained 
linear fit of the Wilson and non-perturbatively improved clover data for the four smaller lattice 
spacings. The value in the chiral limit is obtained from a linear fit (see text for details). 
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